Glassiness Vs. Order in Densely Frustrated Josephson Arrays 
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We carry out extensive Monte Carlo simulations of the Coulomb gas dual to the uniformly frustrated 
two dimensional XY model, for a sequence of frustrations / converging to the irrational (3 — VE)/2. 
We find in these systems a sharp first order equilibrium phase transition to an ordered vortex 
structure at a T c which varies only slightly with /. This ordered vortex structure remains in general 
phase incoherent until a lower vortex pinning transition T p (f) that varies with /. We argue that 
the glassy behaviors reported for this model in earlier simulations are dynamic effects. 
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The glass transition to a frozen disordered state re- 
mains one of the oldest unresolved problems of condensed 
matter physics. While much progress had been made in 
the area of "spin glasses" , in which the glassy state is a 
consequence of intrinsic random frustration in the Hamil- 
tonian, the problem of "structural glasses" , which pos- 
sess no intrinsic random disorder, remains poorly under- 
stood |l],|2| . It is therefore desirable to search for glasslike 
transitions in simple intrinsically disorder-free statistical 
models. One such candidate system is the two dimen- 
sional (2D) uniformly frustrated XY model, which mod- 
els a periodic array of Josephson junctions in a trans- 
verse applied magnetic field [HJJ]. Varying the frustra- 
tion parameter / (magnetic field) of this model leads to 
complex commensurability effects between the underly- 
ing discrete grid and the vortex lattice that forms in re- 
sponse to the frustration Q. Some years ago, Halsey || 
presented numerical evidence that, in the limit of an irra- 
tional /* = (3 — y/5)/2, this model displays a finite tem- 
perature glass transition T g to a superconducting frozen 
disordered vortex state fl. Experiments on supercon- 
ducting wire networks with Halsey's irrational /* have 
found evidence for a finite T g from the scaling of current- 
voltage (IV) characteristics However simulations by 
Granato ||, using resistively shunted junction dynamics, 
found an IV scaling consistent with T g — 0. Recently, 
Kim and Lee have re-investigated this problem using 
Langevin simulations. They find that near Halsey's T g , 
the system's dynamics resembles the primary relaxation 
of supercooled liquids rather than that of a spin glass. 

In view of the above conflicting results, it is important 
to establish the true equilibrium behavior of this model. 
We therefore re-investigate Halsey's problem by carrying 
out Monte Carlo (MC) simulations of the 2D Coulomb 
gas which is dual to the uniformly frustrated XY models. 
Working with vortex variables directly allows us greater 
control over the most relevant slow variables involved in 
equilibration, as compared to using the phase variables 
of the original XY model Following Halsey, 

we consider the frustrations / = 5/13, 8/21, 13/34, and 
21/55, which are the first few members of a Fibonacci se- 



quence of rational approximants which converges to the 
irrational /* = (3 — VE)/2. We find that the low tem- 
perature state that is reached upon slow cooling is highly 
sensitive to both the dynamics used as well as the sys- 
tem length L. The true ground states for such / = p/q 
are much more complex than previously believed, even 
for relatively small values of q. We find that when next- 
nearest-neighbor vortex hops are included in the dynam- 
ics, all cases show clear evidence for a sharp first order 
equilibrium phase transition T c near Halsey's T g , to an 
ordered vortex structure consisting of completely filled, 
completely empty, and partially filled diagonals. Below 
T c , vortices in the partially filled diagonals can remain 
mobile, destroying phase coherence in the direction trans- 
verse to the diagonal. These vortices pin to the grid, 
leading to phase coherence in all directions, only at a 
lower T p {f) that varies with /. We therefore conclude 
that the "glass transition" observed by Halsey, and the 
supercooling observed by Kim and Lee, is a consequence 
of energy barriers in their particular dynamics inhibiting 
what is a true first order equilibrium phase transition to 
a non-glassy ordered state. 

The uniformly frustrated 2D XY model, within the 
Villain approximation [jy], can be mapped [Q onto the 
following Hamiltonian for a one component Coulomb gas 
on a neutralizing background, 



WcgM = \ J> 4 - f)G'(r t - r,)(n, - /) 



(1) 



The sum is over all sites of a 2D periodic square L x L 
grid, m is the integer charge on site i, representing the 
vorticity of the XY phase angle. The frustration param- 
eter /, acting as a background charge density, represents 
the density of flux quanta in the applied magnetic field. 
The interaction is G'(r) = G(r)-G(0), where G(r) is the 
lattice Coulomb potential in 2D, with periodic boundary 
conditions. Charge neutrality, JV = L 2 f , is imposed. 
See Ref. Jl3| for further details. 

The elementary move of our MC procedure is the in- 
sertion of a randomly positioned vortex-antivortex pair, 
Arii = +1, Arij = — 1, which is then cither accepted or 
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rejected by a standard Metropolis algorithm. When we 
restrict i and j to nearest-neighbor sites, we find glassy 
results qualitatively similar to Halsey's. When we allow 
i and j to include next-nearest-neighbor sites as well, we 
find equilibration at low temperatures to be dramatically 
improved. Our results below are for this latter dynam- 
ics. Simulations were carried out cooling from an initial 
random configuration. At each temperature 20,000 ini- 
tial MC passes are discarded for equilibration, with an 
additional 1,280,000 MC passes for computing averages. 
One MC pass refers to 1? elementary moves. 

In Fig. |l] we show intensity plots for the average vor- 
ticity at each site, at T = 0.02 < T c ~ 0.03. Black 
squares are sites with (n*) ~ 0; white squares are sites 
with (rii) ~ 1; gray squares are sites with an average 
vortex occupation of < (rij) < 1. We find an or- 
dered sequence of completely filled, completely empty, 
and partially filled diagonals, clearly different from the 
disordered structures found by Halsey. For Figs. Ek,c,d, 
we find translational invariance along the diagonals, ex- 
cept for occasional defects. For / = 8/21 (Fig. 0b) the 
partially filled diagonals have a pinned vacancy on every 
third site. 



(a) f = 5/13, L = 52 



(b) f = 8/21, L = 42 




FIG. 1. Intensity plot of ordered vortex states at 
T = 0.02 < T c . White squares contain vortices, black squares 
contain no vortices, and gray squares have an average vortex 
occupation of < {rn) < 1. 

We do not know that the states shown in Fig.[l] repre- 
sent the true ordered states in the thermodynamic limit. 
For / = 5/13, for example, using L — 26 and L = 52 
resulted in a differing sequence for the filled, empty, and 
partially filled diagonals. For / = 8/21, the vortices in 
the partially filled diagonals occupy exactly 2/3 of the 
sites in these diagonals. We may speculate that in the 
true ground state these vortices will from a periodic lat- 
tice with the same structure as the / = 2/3 ground state 



Pj . Such a structure can only be made perfectly periodic, 
and commensurate with the periodicity of the diagonals, 
when L is an integer multiple of 84. The structure shown 
in Fig.[|b, with L = 42, consists of such an / = 2/3-like 
arrangement, however with a domain wall introduced by 
our choice of a too small value of L. From such consid- 
erations we conclude that the true ground state, for all 
but the simplest of / = p/q, may involve rather subtle 
and a priori unknown commensurability requirements; 
its square unit cell will be of length mq where m may 
well be an integer several times larger than previously 
believed El. 
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(a) f=5/13,L=52,T c =0.0325 
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FIG. 2. Bimodal energy distribution V(E) at T c . The 
first order transition temperature T c ~ 0.03 is fixed by the 
condition that the peaks subtend equal areas. 

For each case we find that the transition to the ordered 
state is sharp and seemingly first order. To demonstrate 
this, we compute the histogram of energies V(E) encoun- 
tered at each value of temperature in the simulation. In 
both high and low temperature regions, this histogram is 
unimodal. However in a narrow temperature range when 
the ordered state first appears, the distribution becomes 
bimodal corresponding to the two coexisting states at a 
first order transition. Using an extrapolation technique 
p"5| we determine T c as the temperature for which the 
two peaks of the bimodal histogram subtend equal area 
p6[ . Our results are shown in Fig.||. For Figs. |]a,b,d, 
the two peaks are well separated; the less clear case of 
/= 13/34 is perhaps a reflection of the larger concentra- 
tion of defects seen in Fig.|l|c, or the possibility that our 
size L = 68 still gives too poor an approximation of the 
true ground state. In the mapping from the XY model to 
the Coulomb gas the temperature has been rescaled (ljj 
so that T XY = 2ttT cg . Our value T C CG ~ 0.03 is thus 
reasonably close to Halsey's value of T* Y ~ 0.25. 

Next we consider the superconducting phase coherence 
of the Josephson array. In the original uniformly frus- 
trated XY model of phase angles the issue of phase 
coherence may be addressed by considering the depen- 
dence of the total free energy T on the net phase angle 
twist A M that is applied across the sample as a boundary 
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condition, i.e. 8(ri+Lfi) = #(ri) + A M . If J-[A X , A y ] is in- 
dependent of the A^, then phase coherence is lost. Doing 
the duality transformation to the Coulomb gas carefully 
p2| , p7|Jl8| , one finds that such a fixed twist boundary con- 
dition results in an additional term to the Coulomb gas 
Hamiltonian of Eq. (pi), 



V^-A y)+ V { ^ 



SH[p;A x ,A y ] 
where 

P = ^2 VlUi 

i 

is the total "dipole moment" of the vortices, and 



(2) 



(3) 



V(cf>) 



-Tin 



E 



6 4ttT< 



-27trn) 



(4) 



is the Villain function |TT| with coupling 2nT. The par- 
tition function for the system with a fixed twist A p is 
then Z[A x ,A y ] = Z CG (e~ Snlp ' A *' Ay]/T ) ■ Here Z CG is 
the partition function for the ensemble defined by Hcg 
of Eq. (0) alone, and the average is with respect to this en- 
semble (7icc can be considered as the ensemble in which 
A M is averaged over, and so one has "fluctuating twist 
boundary conditions" in the XY model jn]]). The to- 
tal free energy is then F[A x ,A y ] = F C G + SF[A X , A y ], 
where J^cg = — TlnZcc an( l 

^[A*, A,] =-Tln £ V(p x , Py )e- sn ^ A ^/ T (5) 

Px,Py 

where V(p x ,p y ) is the histogram of total dipole moments 
p at a given temperature, found in the simulation using 
Hcg of Eq. By storing this 2D histogram, we can 
therefore deduce the dependence of the free energy on 
all values of applied twist A M . 

In Fig.|| we show intensity plots of 8T[A Xl A y ] for 
A ;t e (— 7T, 7r) at T = 0.02 < T c , corresponding to the 
real space plots of Fig. [l]. For / = 8/21, Fig.||b, we 
see a rotationally symmetric parabolic minimum indicat- 
ing that the system is phase coherent in all directions 
(that there are actually two such minima, is a result of 
the thermal motion of the domain wall inserted by our 
choice of L = 42). For the other cases however, we see 
that while 5 J- has a parabolic minimum along the direc- 
tion parallel to the ordered diagonals, it is constant for 
the direction perpendicular to these diagonals. This in- 
dicates that the vortices in the partially filled diagonals 
are free to move along these diagonals and so destroy 
phase coherence transverse to this direction. Noting that 
an applied electric current exerts a force on the vortices 
which is transverse to the direction of the current, we 
would expect the Joscphson array to have a finite lin- 
ear resistivity for all cases except when the current is 
applied exactly parallel to the partially filled diagonals. 



Thus the structural transition at T c is not in general the 
superconducting transition of the array. 



(a) f = 5/13, L = 52 



(b) f = 8/21, L = 42 




FIG. 3. Intensity plot of total free energy SJ-[A X , A y ] at 
T = 0.02 < T c . Black (white) denotes the functional minima 
(maxima) . 

The helicity modulus tensor T M1/ of the 2D XY model 
is defined by M 



d 2 T 
dAndA„ 



(6) 



where the derivatives are evaluated at the value A M o that 
minimizes T . From Fig. ^ we expect that T M „ is diago- 
nal in a basis that is aligned with the grid diagonal direc- 
tions. We therefore denote by Tj_ and T|| the eigenvalues 
of in the directions perpendicular and parallel to the 
ordered diagonals respectively. In Fig.^ we plot and 
T|| as functions of T, for the same values of / and L as 
in Figs.[l] and ||. As expected from Fig.||| we see that for 
all cases T|| increases from zero as T decreases below T c . 
However, except for the case / = 8/21, remains zero 
below T c , becoming non-zero only at a lower temperature 
T p (f) when the vortices in the partially filled diagonals 
pin to the grid. Similar behavior, of mobile vortex "de- 
fects" in an otherwise ordered vortex structure, has been 
observed previously [^0| in simulations of the / = 5/11 
model. If the vortices in the partially filled diagonals re- 
main mutually correlated, the region T p (f) < T < T c 
can be described as a "floating smectic" phase, as first 
postulated by Ostlund |^l| . If however the correlations 
are short ranged, as in a liquid, one might imagine that 
vortex hopping between the different partially filled di- 
agonals may also be possible. In this case, our result 
that T|| > for T p (f) < T < T c might be a reflection of 
the very high energy barrier for such inter-diagonal hops, 
rather than a true phase coherence effect. 
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FIG. 4. Helicity modulus eigenvalues Tm and Ti vs. T. 
T|| is non-zero below T c ~ 0.03, however Ti is non-zero only 
below T p (f) < T c . 

We conclude that the sequence of rational / that ap- 
proach the irrational /* = (3 — Vo)/2 all undergo a first 
order equilibrium phase transition to an ordered vortex 
structure at a T c ~ 0.03. The exact sequence of the 
filled, empty, and partially filled diagonals in this ordered 
structure remains in general unknown for the true ground 
state in the thermodynamic limit, however in the cases 
when we varied L — mq for fixed / = p/q, we found 
that T c remained approximately 0.03. The transition to 
the true superconducting state, with phase coherence in 
all directions, occurs in general at a lower T p (f), which 
can show considerable variation with the frustration /. 
While this is in qualitative agreement with arguments 
by Teitel and Jayaprakash which suggested that the 
superconducting transition temperature would be a very 
discontinuous function of /, we as yet can discern no sys- 
tematic dependence on / = p/q nor can we be certain 
that the values of T p (f) obtained here will not vary if 
one increases the system size L. 

Our equilibrium transition at T c to an ordered vor- 
tex structure was only obtained when we included next- 
nearest-neighbor hops in our vortex dynamics. When 
moves were restricted to nearest-neighbor hops only, our 
simulations fell out of equilibrium into a frozen disor- 
dered state below T ~ 0.033. Fig. |l| suggests why this is 
so. Next-nearest-neighbor hops allow vortices to travel 
directly up and down the partially filled diagonals. To 
make such a move using only nearest-neighbor hops, one 
must first hop to a neighboring diagonal, putting three 
vortices mutually adjacent. We find for the energy bar- 
rier of such moves AE ~ 0.23 — 0.35, depending on /, 
and so these moves tend to freeze out by the temperature 
T c <C AE. Indeed, when restricting to nearest-neighbor 
hops only, we found glassy behavior even for the simple 
frustration / = 3/8 at large L, unless very slow and care- 
ful cooling was used. Such glassy behavior is therefore 
more a reflection of the frustration being non-trivially 
dense rather than specifically irrational. 

Our results indicate that the "glass transition" ob- 



served by Halsey |^] is an artifact of his choice of dy- 
namics. Our observation of an equilibrium first order 
transition strengthens the analogy to structural glasses 
and gives a natural explanation for the supercooled re- 
laxation observed by Kim and Lee, whose simulations 
were carried out for the parameters / = 13/34, L = 34. 
Further work is required to investigate whether such a 
supercooled state can have a well defined finite temper- 
ature glass transition below T c , or whether, as suggested 
by Granato, this glass transition is at T = 0. 
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